Effect of Layer-Stacking on the Electronic Structure 

of Graphene Nanoribbons 



o 



Neerav Kharche,*'"'''" Yu Zhou,* Kevin P. O'Brien,^ Swastik Kar,§ and Saroj K. 

Nayak*'* 

Computational Center for Nanotechnology Innovations, Rensselaer Polytechnic Institute, Troy, 

NY 12180, USA, Department of Physics, Applied Physics and Astronomy, Rensselaer Polytechnic 

Institute, Troy, NY 12180, USA, Intel Corporation Components Research, Hillsboro, OR 97124, 

USA, and Department of Physics, Northeastern University, Boston, MA 02115, USA 



E-mail: [kharcn@rpi.edut|nay aks@rpi. edu 



^ Abstract 



o 

CZ3 



a 



I 

o 
o 



> 

en 
rn 

O 



> 

X 
J3 



The evolution of electronic structure of graphene 
nanoribbons (GNRs) as a function of the number 
of layers stacked together is investigated using ab 
initio density functional theory (DFT) including 
interlay er van der Waals interactions. Multilayer 
armchair GNRs (AGNRs), similar to single-layer 
AGNRs, exhibit three classes of band gaps de- 
pending on their width. In zigzag GNRs (ZGNRs), 
the geometry relaxation resulting from interlayer 
interactions plays a crucial role in determining 
the magnetic polarization and the band structure. 
The antiferromagnetic (AF) interlayer coupling is 
more stable compared to the ferromagnetic (FM) 
interlayer coupling. ZGNRs with the AF in-layer 
and AF interlayer coupling have a finite band gap 
while ZGNRs with the FM in-layer and AF inter- 
layer coupling do not have a band gap. The ground 
state of the bi-layer ZGNR is non-magnetic with 



*To whom correspondence should be addressed 

^Computational Center for Nanotechnology Innovations, 
Rensselaer Polytechnic Institute, Troy, NY 12180, USA 

^Department of Physics, Applied Physics and Astron- 
omy, Rensselaer Polytechnic Institute, Troy, NY 12180, 
USA 

'"Intel Corporation Components Research, Hillsboro, OR 
97124, USA 

^Department of Physics, Northeastern University, 
Boston, MA 02115, USA 

II Department of Physics, Applied Physics and Astron- 
omy, Rensselaer Polytechnic Institute, Troy, NY 12180, 
USA 



a small but finite band gap. The magnetic order- 
ing is less stable in multilayer ZGNRs compared 
to single-layer ZGNRs. The quasipartcle GW 
corrections are smaller for bilayer GNRs com- 
pared to single-layer GNRs because of the reduced 
Coulomb effects in bilayer GNRs compared to 
single-layer GNRs. 
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Introduction 

Graphene has attracted enormous attention due to 
its extraordinary electronic, optical, thermal and 
mechanical properties, and immense potentials for 
nanoelectronic applications .^^ xhis single-atom 
thick, 5;?^ -hybridized allotrope of carbon with a 
perfectly 2D confinement of its electronic states 
is a zero-gap semi-metal, exhibiting a linear dis- 
persion relation E{k) = hvfk near the meeting 
point of its conical valence and conduction bands 
(the Dirac point) .'^ Graphene can conduct much 
higher current densities than currently used Cu 
interconnects in the integrated circuits (ICs).'^ 



Graphene also exhibits much higher carrier mo- 
bilities compared to the conventional field-effect- 
transistor (FET) channel materials such as Si and 
III- Vs. ^ These unique electronic properties make 
graphene one of the most promising candidate ma- 
terials for both transistors and interconnects in fu- 
ture ICsP 

For many practical applications, such as in 
nanoelectronic devices, graphene needs to be 
patterened into the so-called graphene nanorib- 
bons (GNRs). Typically graphene is patterned 
by selectively removing the material by physi- 
cal etching techniques. Graphene can also be 
patterned by selective chemica l func tionalization 
with Hydrogen^I^nior Fluorine,^^^^ which results 
in graphene nanostructures embedded in function- 
alized graphene: typically a wide bandgap insula- 
tor. The electronic properties of GNRs are very 
sensitive to their width and edge geometry. ^^^^ 
The dependence of edge geometry on the elec- 
tronic structure of GNRs has been mainly inves- 
tigated using theoretical approaces such as the 
tight-binding model and density functional theory 
(DFT). Theoretical studies indicate that the band 
gap of single layer armchair GN Rs (A GNRs) is 
extremely sensitive to their width.ll^n^] Antiferro- 
magnetic ordering at the edges of zigzag GNRs 
(ZGNRs) opens up a band gap while ferromagnet- 
ically ordered ZGNRs do not have a band gap. ^^ ^^ 
The band gap in both AGNRs and ZGNRs can be 
controlled by an external electric field.'^ Exper- 
imental observations of these theoretical predic- 
tions remain ilusive due to the challenging task 
of atomically precise control of the GNR edges. 
The availa ble state-of-the-art transport measure- 
ments'ISnS] are performed on GNRs with width on 
the order of 10 nm, however their edges are rough 
making direct quantitative comparison with theory 
difficult. The experimental techniques, however, 
keep on improving at a rapid pace. For example, 
the magnetic or dering at the edges predicted by 
DFT calculations'^SUIl has been recently observed 
in GNRs with ultrasmooth edges. ^*^ 

Layer stacking provides an additional handle 
to tune the electronic properties of graphene and 
GNRs through interlay er interactions. The ef- 
fect of layer stacking in graphene has been ex- 
tensively studied theoretically .EDESl Recent exper- 
iments have demonstrated tuning of electronic 



properties of graphene by layer stacking. For ex- 
ample, (i) a band gap can be opened up in a bi- 
layer graphene by applying an electric field, ^^ (ii) 
tri-layer graphene behaves like a semimetal in the 
presence of an electric field,'22l and (iii) multilayer 
graphene shows a peculiar conduc tivity spectra de- 
pending on the number of layers .'^^'^ There have 
been several theoretical studies on the electronic 
structure of bilayer GNRs,^'' ^"^ however, similar 
studies beyond bi-layer stacking are lacking except 
Ref.,'23 where the effects of geometry relaxation 
were not included. 

In this paper, we report an investigation of the ef- 
fect of layer-stacking on the electronic structure of 
armchair and zigzag GNRs using ab initio density 
functional theory (DFT). Interlay er van der Waals 
interactions are included to accurately model the 
effects of geometry relaxation on the electronic 
structure. The width and thickness dependence 
of the electronic structure of multilayer AGNRs is 
discussed first, then the energetics of geometry re- 
laxation and magnetic ordering and their effects on 
the electronic structure of multilayer ZGNRs are 
discussed. 

It should be pointed out here that DFT based 
calculations underestimate the band gap of GNRs 
and perturbative correction schemes such as the 
GW method should be used to obtain accurate 
estimates. ^^ GW calculations are computationally 
much more expensive compared to DFT. Since the 
primary goal here is to study the evolution of elec- 
tronic structure as a function of the number of 
layers rather than quantitatively estimate the band 
gaps of stacked GNRs, we have used the compu- 
tationally less expensive DFT scheme for most of 
the calculations. The more accurate GW corrected 
band gaps are calculated for single and bi-layer 
GNRs to provide an estimate of the quasiparticle 
corrections in these nanostructures. 

The atomistic schematics of simulated bi-layer 
AGNRs and ZGNRs with AB Bemal stacking are 
shown in [[figure] [1][]1| ^ a) and |[figure][l][]l| ^b) 
respectively. For three and more layers ABA and 
ABAB... stackings are used. An interlayer dis- 
tance of 3.35 A (same as graphite) is used as a 
starting point in the geometry optimization. Edges 
of GNRs are passivated with hydrogen. Hence- 
forth 1-8 layers thick GNRs are referred as few- 
layer GNRs while multilayer GNR refers to a 
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Figure 1: Schematics of simulated structures. AB 
Bemal stacking in (a) AGNR-7 and (b) ZGNR-8. 
Grey and black colored atoms belong to top and 
bottom layers respectively. The unit cells are de- 
picted by the dashed rectangles. 

GNR obtained by applying periodic boundary con- 
dition to a bilayer GNR in the thickness direction. 

Results and discussion 

The electronic structure of AGNRs is discussed 
first. The band structure of a single-layer and a 
tri-layer AGNR-7 is shown in 
and 



[figure][2][]2;a) 



[figure] [2] [] 2 b) respectively. As the three 



AGNR-7 layers are moved closer to form a tri- 
layer AGNR-7, the 3-fold band degeneracy is 
lifted due to the interlayer coupling and the band 
gap is reduced. Thus the degenerate bands in each 
of the three single-layer AGNR-7 results in three 
non-degenerate bands in a tri-layer AGNR-7. All 
few-layer AGNRs show a similar band structure 



evolution. As shown in [figure][2][]2 e), the band 
gap of AGNR-7 gradually decreases as the num- 
ber of layers increases and shows a tendency to 
approach the multilayer AGNR-7 limit. The com- 
parison of the band structures of tri-layer AGNR- 



7 ( [figure][2][]2 ;b)), AGNR-8 ( [figure][2][]2 c)) 



and AGNR-9 ( [[figure] [2] []2t d)) illustrates the fact 
that similar to single-layer AGNRs tri-layer AG- 
NRs show three classes {i.e., N = 3p, 3p+l, and 
3p + 2, where A^ denotes the number of carbon 
chains along the width and p is an integer) of 
band structures depending on their width. Similar 
to AGNR-7, AGNR-8 and AGNR-9 show conver- 



gence of the band gap to the multilayer limit ( [fig- 
ure][2][]2) ^e)). This classification was recently re- 
ported in Ref.'23 for upto 4-layer thick AGNRs. 
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Figure 2: Band structure of single- and tri-layer 
AGNRs belonging to three classes 3p, 3p + I, 
and 3p + 2 depending on their width, (a) Single- 
layer AGNR-7, (b) tri-layer AGNR-7, (c) tri-layer 
AGNR-8, and (d) tri-layer AGNR-9. (e) DFT band 
gaps (open symbols) of AGNRs as a function of 
the number of layers. The horizontal lines show 
the DFT band gaps of multilayer AGNRs obtained 
by applying periodic boundary conditions in the 
thickness direction. The solid symbols are the GW 
corrected band gaps of single- and bilayer AGNRs. 



[figure] [2] []2[ e) indicates that this classification 
persists beyond 4 layers. 
The GW corrected band gaps of single- and bi- 



layer AGNRs are also shown in [figure][2][]2 e). 
Large quasiparticle corrections on the order of 
1 — 2 eV are attributed to enhanced Coulomb 
effects. ^^ Compared to bulk, GNRs are under 
confinement and have greatly reduced screening, 
which enhance the Coulomb effects. The quasi- 
particle corrections for bilayer AGNRs are smaller 
compared to single-layer AGNRs. This is be- 
cause each layer screens the Coulomb interactions 

in another layer reducing the overall Coulomb ef- 
fects.SSSQ] 

ZGNRs show an edge magnetism, which 
can have multiple configurations depending on 
relative in-laye r and interlayer spin polariza- 
tions .1IM 18|31}|351 jj^g (^^Q possible spin polariza- 
tions in a single-layer ZGNR are (i) ferromagnetic 
(FM) in-layer and (ii) antiferromagnetic (AF) in- 



layer while multilayer ZGNRs can have several 
possible spin polarizations. The four spin polar- 
izations in multilayer ZGNRs investigated here 
are (i) FM-FM: FM in-layer and FM interlayer, 
(ii) AF-FM: AF in-layer and FM interlayer, (iii) 
FM-AF: FM in-layer and AF interlayer, and (iv) 
AF-AF: AF in-layer and AF interlayer. Multi- 
layer ZGNRs initialized in configurations (iii) and 
(iv) stay in those configurations while multilayer 
ZGNRs initialized in configurations (i) and (ii) 
converge to configuration (iv), which indicates 
that the FM interlayer coupling is not as stable as 
the AF interlayer coupling. Ref. ^^' also reported 
that the AF interlayer coupling has lower energy 
compared to the FM interlayer coupling. There- 
fore, AF interlayer coupling is used in all ZGNR 
electronic structure calculations presented in this 
article. 

The relaxation energy (i.e. difference between 
total energy of relaxed and unrelaxed) ZGNRs as a 



function of the number of layers is plotted in [fig 



ure][3][]3| ^a). Both non-magnetic (NM) and mag- 
netically ordered ZGNRs show similar trends in 
the relaxation energy as a function of the num- 
ber of layers. Due to the interlayer van der Waals 
attraction, the top and bottom layers show a con- 
cave curvature while the other layers remain more 
or less flat. Similar curvatures were reported in 
Ref.^^ for the bi-layer ZGNR. This concave cur- 
vature is the main cause of the reduction in total 
energy of the relaxed few-layer ZGNRs. The re- 
laxation energy of the single-layer ZGNR is negli- 
gible because of the absence of curvature. 



[figure][3][]3;b) shows the magnetic stabiliza- 
tion energy (AEm = E^m — Em) of unrelaxed 
ZGNRs resulting only from the FM-AF and AF- 
AF orderings at the edges. ZGNRs with AF-AF 
and FM-AF orderings both show similar trends 
as a function of the number of layers such that 
ZGNRs with AF in-layer configuration are lower 
in energy compared to ZGNRs with FM in-layer 
coupling. AEm is much smaller in ZGNRs with 
two or more layers compared to the single-layer 
ZGNR. This implies that the interlayer interations 
weaken the strength of magnetic ordering in each 
layer of the few-layer ZGNRs.^H The interlayer 
interactions become stronger in relaxed ZGNRs 
because of the concave curvature in the top and 
bottom layers. This further weakens the strength 
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Figure 3: Energetics of ZGNR- 8. Difference 
between total energy of (a) relaxed and unre- 
laxed ZGNRs, (b) magnetically ordered and non- 
magnetic unrelaxed ZGNRs, and (c) magnetically 
ordered and non-magnetic relaxed ZGNRs as a 
function of the number of layers. 

of in-layer magnetic coupling, which results in 
lower energy difference between magnetic and 
non-magnetic states of relaxed ZGNRs compared 



to unrelaxed ZGNRs [[figure] [3] []3i :b,c). 

The bi-layer ZGNRs initialized in both the FM- 
AF and AF-AF configurations converge to a non- 
magnetic ground state when the geometry relax- 
ation is allowed. Therefore in the bi-layer ZGNR, 
there is no reducti on in the total energy due to the 
magnetic ordering [figure] [3] [] 3 c). Ref.,^^ which 
included van der Waals interactions also reported 
that the ground state of bi-layer ZGNRs is non- 
magnetic. In the calculations, where the geometry 
relaxation is not included, the bi-layer ZGNRs ini- 
tialized in the FM-AF and AF-AF configurations 
stay in those configurations and do not converge 
to a non-magnetic ground state [figure] [3] []3[ b). 



The calculated AEm of relaxed ZGNRs ([fig 



ure][3][]3 c)) is comparable to the thermal energy 
at room temperature (fcgT ^ 25 meV), which in- 
dicates that magnetic ordering is unstable at room 
temperature. The layer stacking destabilizes edge 
magnetism by reducing AEm- The roughness and 
reconstruction at the ZGNR edges further desta- 
bilize the magnetic ordering even at low temper- 
atures. Thus the magnetic ordering is expected 
to manifest at low temperatures in single-layer 
ZGNRs with atomically smooth edges. Indeed the 
edge magnetism was observed recently at low tem- 
peratures in GNRs with atomically smooth edges 
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Figure 4: DFTband structures of 1-4 layer ZGNR- 
8 in non-magnetic (a-d), antiferromagnetic (AF) 
in-layer and AF interlay er (e-h), and ferromag- 
netic (FM) in-layer and AF interlayer (i-1) config- 
urations. Single-layer ZGNRs do not have AF in- 
terlayer couplings. Inset in (h) shows a zoomed in 
view of the band structure near the band gap. The 
dotted lines in (e) and (f) are the GW corrected 
band structures. The up and down spin states in 
(e-h) are degenerate while they are non-degenerate 
in (i-1). 

slightly offset from an ideal zig-zag edge-I^SI 

The band structures of 1-4 layer thick relaxed 
ZGNRs in the non-magnetic, AF-AF, and FM-AF 



configurations are shown in [figure] [4] []4 Non- 
magnetic ZGNRs have nearly fiat bands near the 
Fermi level (Ep) leading to high density of states 
(DOS) near Ep. The bands near Ep are mainly 
composed of the edge states. The high DOS near 
Ep gives rise to magnetic instability and the edge 
states become spin-polarized. ^^ The up and down 
spin states in ZGNRs with AF-AF ordering at 
the edges are degenerate while they are split in 
ZGNRs with FM-AF ordering. The AF in-layer 
coupling at the edges opens up a band gap while 
FM in-layer coupling does not. The band struc- 
tures of 5-8 layer thick relaxed ZGNRs are not 
shown here but are qualitatively similar. The band 



structures of the bi-layer ZGNR in the AF-AF 
( |[figure][4][]4| :f)) and FM-AF ( [[figure] [4] []40 )) 
configurations are identical to non-magnetic band 
structure ( [[figure] [4] [] 4 b)) because the bi-layer 
ZGNRs initialized in those configurations con- 
verge to the non-magnetic ground state after ge- 
ometry relaxation. 

In the presence of the magnetic ordering at the 
edges, the non-magnetic bands in ZGNRs are only 
slightly perturbed except the bands near Ep. Al- 
though the band structures of ZGNRs with even 
and odd number of layers look almost similar, 
they are different near Ep and the Brillouin zone 
edge. Non-magnetic ZGNRs with odd number of 



layers have nearly flat bands crossing Ep ( [fig 
ure][4][]4[ ^a,c)). Such bands are not present in 



non-magnetic ZGNRs with even number of lay- 
ers. This property is reminiscent of the differ- 
ent electronic structure of graphene multilayers, 
where Dirac fermions with a linear dispersion are 
present in graphene with odd number of layers 
while only normal fermions with a parabolic dis- 
persion are present in graphene with even number 
of layers. ^^ Upon magnetic ordering at the ZGNR 
edges, the splitting of the bands crossing Ep is 
higher compared to the other flat bands slightly 
away from Ep. The bandstructures of ZGNRs 
with 4 or more layers are in close agreement with 
each other irrespective of their magnetic order- 



ing ([figure] [4] []4^d,h,l)). This is also reflected in 



their total energies in [figure][3][]3 c), where the 
energy difference between non-magnetic, AF-AF, 
and FM-AF orderings become very small as the 
number of layers increases beyond 3. 

The band gap of ZGNRs as a function of the 
number of layers is shown in [[figure] [5] []5[ As 



discussed earlier, only ZGNRs with AF-AF or- 
dering and a non-magnetic bilayer ZGNR have a 
band gap. Relaxed ZGNRs have a smaller band 
gap compared to unrelaxed ZGNRs. The GW cor- 
rected band structures of a relaxed single-layer 
ZGNR with AF in-layer ordering and relaxed bi- 
layer ZGNR, which has a non-magnetic ground 
state are shown in [figure] [4] []4| ^e,f). The GW 
corrected band gaps of single- and bilayer ZGNRs 
also shown in [[figure] [5] []5[ The quasi- 



are 



particle corrections for unrelaxed bilayer AF-AF 
ZGNR are smaller compared to unrelaxed single- 
layer ZGNR with AF in-layer ordering because 
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Figure 5: DFT band gaps (open symbols) of re- 
laxed and unrelaxed ZGNRs in AF-AF configura- 
tion as a function of the number of layers. The 
band gap of relaxed bi-layer ZGNR, which has a 
non-magnetic ground state, is depicted by an open 
square. The corresponding GW corrected band 
gaps of single- and bilayer ZGNRs are depicted 
by the solid symbols. 

each layer screens the Coulomb interactions in 
another layer reducing the overall Coulomb ef- 

a non-magnetic relaxed bilayer ZGNR are much 
smaller compared to unrelaxed bilayer ZGNR with 
AF-AF ordering. 

Since the electronic structure of ZGNRs is 
mainly determined by the edge states ™ ZGNRs 
with different widths are expected to show similar 
electronic structure evolution with the number of 
layers as ZGNR-8 discussed here. Ultra-narrow 
ZGNRs may show different electronic structure 
evolution due to strong interedge interactions. 

Conclusion 

To summarize, the evolution of electronic proper- 
ties of GNRs as a function of the number of layers 
stacked together is studied using DFT including 
van der Waals interactions. Multilayer AGNRs, 
similar to single-layer AGNRs, are found to ex- 
hibit three classes of band gaps depending on their 
width. In ZGNRs, the AF interlayer coupling is 
more stable compared to the FM interlayer cou- 
pling. ZGNRs with the FM in-layer and AF in- 



terlayer coupling do not have a band gap while 
ZGNRs with the AF in-layer and AF interlayer 
coupling have a finite band gap, which decreases 
with increasing the number of layers. The mag- 
netic stabilization energy decreases as the num- 
ber of layers increases indicating that the magnetic 
ordering is less stable in multilayer ZGNRs com- 
pared to single-layer ZGNRs. The ground state of 
the bi-layer ZGNR is found to be non-magnetic 
with a small but finite band gap. The DFT calcu- 
lations, which do not include geometry relaxation 
can not predict the non-magnetic ground state of a 
bi-layer ZGNR and overestimate the band gap of 
multilayer ZGNRs in AF-AF configuration. The 
GW calculations on single- and bilayer GNRs in- 
dicate that the quasiparticle band gap corrections 
decrease with increasing number of layers due to 
the reduction in Coulomb effects. 



Methods 

The electronic structure calculations are per- 
formed within the framework of first-principles 
density functional theory as implemented in the 
Vienn a a b initio simulation package (VASP) 
code. ^^^ The PAW pseudopotentials^^^^ and the 
PBE exchange correlation functional in the gen- 
eralized gradient approximation^^ are used. The 
DFT-D2 method of Grimmel^ as implemented in 
VASP^^ is used to model the van der Waals (vdW) 
interaction between GNR layers. To ensure negli- 
gible interaction between periodic images, a large 
value (10 A) of the vacuum region is used. The 
ID Brillouin zone of few-layer GNRs is sampled 
using 32 uniformly spaced ^-points while the 2D 
Brillouin zone of a multilayer GNR is sampled 
using 1 X 16 X 10 Monkhorst-Pack mesh.HSlFor the 
plane wave expansion of the wavefunction a 400 
eV kinetic energy cut-off is used. The total energy 
and the atomic force are converged to within 10^'^ 
eV and 0.05 eV/A respectively. To obtain the band 
structure of few-layer GNRs a non- self-consistent 
calculation is carried out on 101 uniformly spaced 
^-points in the positive half of the Brillouin zone 
using the converged charge density from the self- 
consistent calculation. 

The GW calculations are performed using the 
ABINIT code. — The norm-conserving pseudopo- 



tentials generated using the Trouiller-Martins 
scheme ^^ implemented in the fhi98PP pseudopo- 
tential program ^^ are used. The PBE parameter- 
ization for the exchange-correlation functional"^ 
is used. To ensure negligible interaction between 
periodic images, a large value (10 A) of the vac- 
uum region is used. The Brillouin zone is sampled 
using 32 uniformly spaced ^-points. For the plane 
wave expansion of the wavefunction, a 12 Ha 
kinetic energy cut-off is used. The DFT band 
structures calculated using VASP and ABINIT are 
virtually identical. The quasiparticle corrections 
are calculated within the GqWo approximation 
and the screening is calculated using the plasmon- 
pole modeLlSSThe Coulomb cut-off technique pro- 
posed by S. I. Beigi et alP^is used to minimize 
the spurious interactions with periodic replicas of 
the system. 
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